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ABSTRACT: 

Fluid flow and heat transport are investigated for the GaAs melt rotating with its cylindrical 
container in a liquid encapsulated Czochraslki (LEC) growth system in which the flow pattern and 
temperature distribution in the melt are influenced by coupling of rotationally-driven forces to the 
thermal convection. The experimentally observed particular shape of the crystallization front was 
predicted numerically. The convective flow in the melt was studied by focusing on the Rayleigh-Benard 
and the baroclinic instabilities. Increasing the crucible rotation rate, the dominant convective cell found 
to pivot by a small angle in the meridional plane around the vorticity vector. This pivoting angle is 
directly related to the Ekman shear angle. 
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I. INTRODUCTION 

Low Prandtl-number ( ) semiconductor melt in a modified Czochralski (Cz) puller [1, 2] can be 
modeled as a two-dimensional Boussinesq fluid which rotates uniformly with its container. In the absence of an 
applied electromagnetic field, hydrodynamic processes relevant to crystal growth are classified into the 
following main groups: natural gravitational type convection, forced convection and the surface-tension driven 
convection. Typically, hydrodynamics and the associate transport phenomena in the system are quite complex 
because of their highly nonlinear and strong coupled interactions [3, 4]. The crystal quality is closely related to 
the crystallization front shape which depends on the heat balance at the vicinity of the front. Associated with a 
curved (convex or concave to melt) phase interface, the radial non -uniformity of the thermal history affects the 
crystal structure due to the formation of facets on the oxides growth interface [5]. In the case of semiconductors 
such as GaAs and InP, the density of the dislocations is a direct consequence of thermal stress in the crystal 
during the process. The maximum stresses which occur at the crystal periphery were found [6, 7] to be related to 
the growth interface morphology. The gull-winged shape GaAs phase boundary, revealed by different model 
calculations [8, 9], is well supported by experimental data [10, 11]. This typical interface shape, concave to melt 
( ) at periphery and convex ( ) to melt in the center, strongly depends on the flow field structure [12]. Therefore, 
the interface destabilization, leading to thermal stresses and high dislocation density in the crystal, is mainly due 
to the convective instabilities in the melt. 

Hydrodynamic stability, concerning the stability and instability of the fluid motions, has been largely 
studied in the last decades [13, 14], Experimental data on the regime diagram of different fluids, each contained 
in a rotating cylindrical annulus and subjected to an imposed radial temperature gradient, were reported first by 
Fein and Pfeffer [15], Their results were summarized in the form of two -parameter diagrams (thermal Rossby 
number vs. Taylor number) for different Pr numbers. For moderate and low Pr-number fluids, namely water ( ) 
and mercury ( ), they have displayed the well-known knee type (or the anvil shape) diagrams which characterize 
the transition from axisymmetric flow to regular baroclinic waves. The anvil-shape regime diagram for 
moderate Pr-number fluid, apart from the baroclinic regular wave, revealed the upper and the lower symmetric 
regimes. While the flow in the lower one is established by diffuse effects, the axisymmetric flow in the upper 
symmetric regime is stable due to stratification. 

Though geometry of Cz configuration is different from the rotating annulus experiments, similarity of the flow 
field structure in between has been paid much attention [16, 17]. Compared to moderate and high Pr-number fluids, studies 
on hydrodynamic stability of low Pr-number Cz-melts (Si, GaAs, InP), initiated by Ristorcelli and Lumley [16] and 
Kakimoto [17], appears to be scare. This is more pronounced in the case of GaAs melt for which the convective flow was 
studied by Polezahaev et al [18], Using an axisymmetric modeling of fluid flow without a swirl, the critical Grashof number 
was estimated to be as high as - for GaAs/Cz system of the aspect ratio [18]. They have shown that the critical number for 
the onset of thermal oscillations in the melt decreases to - when the melt depth increased, namely for . However, it is 
well-known [19] that for a low Pr melt rotating with its container, the critical Rayleigh number increases with the rotation 
rate. The so-called overstability associated with rotation intensity, namely the Taylor number, is an important feature of 
convection with rotation [20] . 
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The main purpose of the present work is to investigate the flow field structure of the GaAs/Cz melt. 
The effect of rotationally driven- and thermocapillary forces on the flow pattern is studied. Both rotating 
Rayleigh-Benard instability [21, 22] and the baroclinic instability [17, 23] are to be considered. The 
mathematical model developed in the present study incorporates transport processes of all domains of the 
system. The crystal is considered to be passive here and serves the purpose of enforcing an isothermal boundary 
condition at the crystallization front. 

II. MODEL DESCRIPTION AND NUMERICAL SCHEME 

2.1. Physical Model and Basic Assumptions 

The idealized geometry of a modified Cz puller, the computational domains and the coordinate system 
adopted in the present analysis are illustrated in figure 1. The system is assumed to be in a pseudo -steady state 
and the axi-symmetric. Continuity, Navier-Stokes and energy equations are solved with the Boussinesq 
approximation. Both the GaAs melt and the encapsulant are incompressible Newtonian fluids with laminar 
flows. Thermophysical properties of the system are constant except for the temperature-dependency of the 
encapsulant dynamic viscosity, p e . The crucible bottom is thermally insulated and its side wall of height (h m +h e ) 
is at a uniform and constant temperature, T c . The physical properties of the system used for simulations are 
given in reference [24]. The thermophysical properties of ambient gas (Argon) are summarized in Table 1. The 
geometrical as well as the process parameters used in the present simulations are given in Table 2. 



Table 1: Physical properties of Gas used for simulations 



Density, p 


486.8295/T kg/m 3 


Thermal conductivity, A g 


0.01+(2.5xl0* 5 ) xT W/m.K 


Specific heat, c Pg 
Dynamic Viscosity, p 


521 J/kg.K 

8.466xl0~ 6 + 5.365x10 s T- 8.682xl0 42 T 2 



Table 2: Geometrical and process parameters used for simulation 



Crucible radius, r c 


0.075 m 


Crystal radius r s 


0.0375 m 


Melt height, h„, 


0.05 m 


Encapsulant height, h e 


0.02 m 


Insulating enclosure height, hE 


0.2 m 


Characteristic length 


L=[h 4 m /(r c -r s >] 1/3 


Crystal height, h s 


0.12 m 


Crystal rotational rate, Q x 


5 rpm 


Crucible rotational rate, -Q c 


5-20 rpm 


Pulling rate, u p 


10 mm/h 


Crucible wall temperature, T c 


1529 K 


Enclosure top temperature, T a 


1261 K 


Insulating enclosure emissivity, £e 


0.8 



In this configuration, the counter rotation of the crystal (r s ;Qj and crucible (r c : - n c ) influences 
the flow field in the melt. In the present study, a set of numerical simulations was performed for a 37.5mm 
radius crystal, 75mm radius crucible and 50mm height melt. Hence, the ratio between the Reynolds numbers 

I I 2 

Re c and Re s equals to 4fl ( /n J by definition of Re = r n / v . Throughout the calculations, O x is assumed 
to be constant and equal to a reasonably [25] small value (5 ipm = 0.524 rad / v) while the crucible rotation rate 
varies from -5 to -20 rpm. 
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Figure 1: Sketch of a modified Czochralski crystal growth setup and the computational domains. The curve PIJ 
represents the M-type experimental interface shape on which the point I is the inflection point. 

2.2. Mathematical Model 

Throughout the present simulations, the Reynolds number of the GaAs melt flow field R e a = « lm L/v, 

does not exceed the computationally obtained [25] Re x and/or Re c values. This may ensure the validity of the 
axisymmetric flow assumption in the present modeling. Therefore, all variables o> = <t>(r,z) in Table 3 are 
independent of the azimuthal 0 angle and the numerical model solves three momentum equations for three 
velocity components , « r , «, and u 0 in two-dimensional space. The variations in the fluid density are ignored 
except insofar as their effect on the gravitational force. 



Table 3: The diffusion coefficient r 0 and source term s$, for variable <I> 



<l> 


G 


s. 


1 


0 


0 


u 

z 


M 


— (3 P / Bz) ■+■ pg /3 (/ — T re/ ) 


u . 




— ( dp / dr ) — jli u r / r 2 + /? u ~ / r 


U 6 


M 


— // w q / r — p u r u q / r 



Based on these assumptions, the generalized transport equation (GTE), governing the fluid flow and heat 
transfer in this arrangement, can be expressed in polar-cylindrical coordinate system as 
8 , . 1 8 . . 0 ( 3 ®') 1 ( 3 ® 3 

(p«-®)+ ( pru r ®)= 1 r 0 I + — I rV (D I + s® , (1) 

dz r dr dz k 8z ) r ( dr ) 

for which the diffusion coefficient r 0 and the source term s ([) for dependent variables <t> are listed in Table 3. 
The first and the last rows in Table 3 are the expressions of continuity and the energy equations, respectively. 
The energy loss due to viscous dissipation is neglected. 7 rf , ( is the reference temperature for the Boussinesq 

approximation, p the reference pressure, P the coefficient of volumetric expansion and g is the gravitational 
acceleration constant. 

In the present cylindrical model, the streamline function \|/ and vorticity co are defined in terms of the 
radial and axial velocities as u r = -(1 / r)( d y/ / dz) , u, = {Hr)(dy/ I dr) and a> = du r / dz - du / dr . 
Physically \|/ is a measure of the volume flux. So the difference between the values of t|/ at two points in the 
meridional plane is proportional to the meridional flux between them. The vorticity is associated with rotational 
motion and in a 2D flow the vorticity vector is orthogonal to the plane of flow. 
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2.3. Boundary conditions 

Corresponding to the boundaries defined in figure 1, the conditions applied to the system are listed in Table 4. 
Altogether, three types of flow boundary conditions were used for the whole set of computations: no slip boundary for the 
crucible, crystal and insulator, no shear stress on symmetry axis and Marangoni boundary condition for both the encapsulant 
free surface and the melt/encapsulant interface. The numerically approximated temperature [24] at the crucible side wall ( 
0 < z < (h m + h e ) ) was employed as a constant. Corresponding to the heat input around 1.5 kW in the case of mixed 
convection, T c =1529 K is suitable for growth of 3" diameter GaAs crystal providing the pull rate ~ 10-15 mm/h in the 
growth process. The crucible bottom is insulated and assumed to be adiabatic. Above the crucible top end, the temperature of 
the insulating enclosure wall, T E is assumed to be variant as given in Table 4. The top enclosing wall is at the same 
temperature as the ambient T a =1261 K. Both the temperatures along the crucible bottom and the enclosure top surface are of 
considerable importance for the control of the GaAs/Cz process, particularly for the shape of CMI. 



OP 0 

AP 0 

AB 

AM 

MC 

DC 

KD 

FK 

OF 

JK 

JG 

GB 

GD 



Table 4: The applied boundary conditions corresponding to Figure 1 . 



Ur 



UP 



5m . 
Or 






d u r 
dz 



0 

0 

0 

8r„, ST m 

8T 8r 

0 






d u r d y e d T e 



0 

0 

0 

r c n c 
r c n c 
r « c 

Su e _ 
dz 

r x^ . 
r O 



Sr 

dT 

dr 



dT 44 

-= h se ( T s - T a ) +s s ct(T * - T ) 



dz 



sg v s a 

dT 
dr 



r h E -h z i 

T 7 =(T c -T„ )| — “I + T„ 



eff 



’nrM 



— a 



- A 



3 T m 
dz 

ar 

dr 

dT 



T c 
T c 
dT 



dT 4 4 

4- A — = oc s cr (T — T ... ) 

dz 

dT p 44 

A „ = cc £ ex {T — T ... ) 



dr 



dr 



= h (T - T )+ £ cr (T - T _ ) 

se v s a ' s v s eff ' 



dz dT dr 



dT 44 

- A e—^= h eg (T e -T a ) + e e cr^Te ~ T ) 
dz JJ 



The effective temperature, T e ff, k which appears in Table 4, is due to radiative heat exchange in the present 
system. According to the Gebhart theory of radiation, the net radiative heat exchange between an exposed 
surface k at temperature T k , and all others in an N-sided enclosure can be calculated as: 

N 

9k = T* - — ^ G ik SiO T, 4 A t = e k a (T 4 - k ) (2) 

i.l 

Where A k is the area of the surface k and 

N 

rjt , ^ X™ 1 A /-• ^4,1/4 

T eff ,k ~( X A i £ i G * T i ) (3) 

Ak£k i = \ 

is the effective temperature to which the surface k loses heat by radiation. The Gebhart factor, 

N 

G ik = F ik s k +£ Fij(l- £j )G jk ( 4 ) 

j = 1 

represents the fraction of out-going flux from surface i that absorbed by surface k. The calculation procedure for 
the view factor matrix F(i,j) is straightforward for meniscus free surface. The above set of linear equations for 
G ik (i, k =1,2,...,N) can be solved by Gaussian elimination. Repeating the calculations for all surfaces by 
changing k, one obtains all the Gebhart factors. For N=5 flat-zones enclosure in the present configuration, the 
number of independent view factors to be determined is only (N/2)(N-3)=5. 
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2.4. Numerical model 

The two-dimensional differential equations (1) with the boundary conditions for the fluid flow and heat 
transport in the Cz configuration were solved numerically by the finite volume method (FVM). The 
computational domain was discretized into a finite number of control volumes over which the governing 
equations were integrated. The resulting system of algebraic equations was solved iteratively until convergence 
was reached. A staggered grid system was employed for different dependent variables. The pressure p, the 
azimuthal velocity u e , and the temperature T were calculated at the grid point in the center of the main control 
volume (MCV); the axial and radial components of velocity «_ and « r , respectively, were staggered in the z 
and r directions to the mid-point of the boundaries of MCV. The SIMPLEC algorithm [26] was used to couple 
velocities and pressure on staggered grids. The grid used in the present calculation is uniform for the 
computational domain of melt, encapsulant, crystal and the ambient gas. It contains 52 x 272 grid points with 
spacing of 1.0 mm resulting in a grid-independent solution. For all variables <f>, the solution was deemed 



convergent when the criterion 




- <t> 



^ to 4 was satisfied where n represents the index of iteration 

max 



number. The accuracy of the numerical code was ascertained by validating the general results of the present 
calculation with the experimental [11, 27] and 2D numerical [6, 9] results obtained for corresponding geometry 
and boundary conditions. Remarkable is the particular shape of CMI which appears in the present simulations as 
obtained in the model analyses of both liquid encapsulated Czochralski and VCz growth of GaAs bulk crystals. 



III. RESULTS AND DISCUSSION 

In this section, the nature of the GaAs/Cz melt convective flow is described. The results obtained from 
the numerical simulations are presented. The influence of rotationally-driven forces on the flow field structure 
and temperature distribution is discussed. 



3.1. Description of the flow field structure 

The convection instabilities in Czochralski melt are among the most challenging issues in bulk crystal 
growth modeling. This is mainly due to the complexity of hydrodynamics and associated heat transport in the 
melt from which a crystal is growing. In the present model, the flow instability can exhibit unique features 
because of coexisting vertical nad horizontal temperature gradients and the differential rotation rates of the 
crystal ( r x ; ft x ) and crucible ( r c ; - n c ). The minus sign here means the counter-clockwise direction of rotation. 
For a liquid rotating with its container, besides the gravity, it is necessary to include two body forces, the 
centrifugal force representing a vector derivative £ 2 C a (Si c at) and the Coriolis force 2 (t > a u) with r and 
u as the position and velocity vectors, respectively. The effects of centrifugal force are usually neglected in 
model calculations, particularly for the crystal growth set-ups, and the analysis has largely been restricted to 
consideration of the interaction between the Coriolis and buoyancy forces. In the GaAs melt, the natural 
convection may readily become unsteady and time -dependent because of its low Pr number. Even with small 
driving temperature differences, a r max = 18 K imposed on the system (h m / r c = 0.67 ), the Grashof number 
Gr m was found to be of order 10 7 . However, compared to its critical number calculated for 2D melt without a 
swirl [18], the flow may be assumed to be axisymmetric in the present model in which the melt is rotating with 
its container. 

The crystal rotation generates an upward flow in the melt central column which encounters the large- 
scale circulation (LSC) caused by the crucible rotation in the same direction as the buoyancy and 
thermocapillary convections. In the configuration with /u < o as the rotation ratio -q c / n t and the radius ratio 
of t] — r x I r c — 0.5 , this can be readily verified that Rc c / Rc x =4/; . With a constant and reasonably small [25] 
n x = 0.524 rad / ,v throughout the calculations, the large-scale flow was found to be the dominant feature of the 

flow pattern. The parameter Ra / Pr .Re 7 is so high (~10), compared to Ra / Pr . Re 7 that the upward central 
flow does not appear in the present simulations. In fact, only for n < 0.417 rad / 4 (~4.0 rpm), the governing 

parameter (Ra / Pr . Re 7 ) > 1 can be held. Hence, the reason for an eventual instability in the system is coupling 
of thermal convection and the forces generated by the crucible rotation. Furthermore, the surface tension to 

buoyancy ratio Ma / Pr .a 1 7 2 for the GaAs melt is large enough (-180 compared to -120 in[16] for Si) to 
modify the structure of the flow. However, as illustrated in figure 2, the effect is strongly controlled by the 

presence of highly viscous (/j e //r m ~ to 3 ) encapsulated layer on the GaAs melt surface. Far from the 
melt/encapsulant boundary, the flow velocity vectors and temperature profiles in the melt are influenced by 
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thermocapillary forces. This effect leads to a more convex to melt interface by changing slightly the direction of 
the velocity vectors in peripheral part of the phase boundary. The flow velocity vectors close to the 
melt/encapsulant interface were found to be so small (see figure 2) and the flow remained in a steady state as 
assumed in the present model. This can be inferred that, for the GaAs melt/encapsulant system, the effective 

critical Ma number should be far beyond the estimated - 6.9 x 10 4 (Pr) 1 36 value [28] for low Pr melts. 




Figure 2: Marangoni effect (black lines) on the velocity field (left) and temperature distribution (right). The 
rotation rates for crystal and crucible are +5 and -5 rpm, respectively. 

For the GaAs melt in a rotating Cz crucible, the Coriolis force largely influences the buoyancy-driven 
flow due to coupling of the radial velocity with the fluid rotation. Though the Coriolis force stabilizes the fluid 
flow due to the prevention of vertical convection, it complicates however the structure of the flow field by 
increasing the curvature the fluid particle-paths. Therefore, hydrodynamics of the GaAs melt in a non-inertial 
frame is characterized by the Taylor number T a and thermal Rossby number Ro T . By definition Ta stands for the 
ratio of Coriolis force to the viscous force, and the Rossby number is the ratio between convective acceleration 
and the Coriolis force. Non-dimensional parameters related to the present model are given in Table 5. 



Table 5. Dimensionless similarity parameters for the GaAs melt convection with the rotation. The characteristic length 
L = 5.503 x 10 ^ 2 m , A T lmx = 18 K , n x = 0.524 rad Is and Q x < n c < 2.094 rad Is. 



Ra = (gfl Iva ) L 3 Ar mx 


1 .574 x 10 6 


Ma =(r r //dz)tA7' nK 


5.928 x 10 4 


Re * = r x n x /v 


1 .511 x 10 3 


R e c = r c 2 Q c / v 


i.i53 x io 4 n c 


2 5 2 2 

Ta =4 Cl c (r c - r x f /v L 


2.744 x 10 7 n 2 


Ro T = sfi A 7-^ L / n 2 c (r c - r x ) 2 


1.292 / n 2 
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(d) Q c = -12.5 rpm 



(e) Q c = -15 rpm 



(f) Q c = -20 rpm 



Figure 3 (a-g): Contours of stream function (left) and temperature (right) with Ra = 1.574 x 10 6 , Ma = 5.928 x 10 4 

, Re x = 1 .511 x 10 3 for different crucible rotation rates. 







Convection interaction in GaAs/LEC... 




Figure 4: Maximum intensity of the convection as a function of the rotation rate. 




Figure 5: Maximum vorticity of the flow as a function of the rotation rate. 

3.2. Simulation of the flow field and temperature distribution 

In the GaAs melt rotating with its container, the flow field and heat transport are influenced by the 
impact of Coriolis force on the thermal convection flow. The flow field and temperature distribution are 

illustrated in figure 3 (a-g) for the melt characterized byRa m = 1.574 xio 6 , Ma = 5.928 x 10 4 , 

3 I I 

Re x = l .511 x to . The crucible rotation varies in the range 0.524 < n c < 2.094 rad is. To describe the multi- 
cell structure of the flow, the line b ( r = r x /2 , z = h m ), crossing approximately the inflection point ( I) on the 
phase boundary, and the line E as the edge of the melt central column (see figure 1) are used. 

As for any rotating stratified fluid in a container [20] the structure of the flow is dominated by a 
centrifugally-driven large-scale circulation (LSC). The roll (LSC) is characterized by i// ( z) I y/ max = l and 
found to be headed towards the peripheral concave to melt part of the phase interface. For all cases, the central 
position of LSC lies approximately on the edge E nearby the point M ; ( r = r c /2 , z = h m / 2 ) in the interior of 

the melt. As shown in figure 3 (a-g), at lower rotation rates |n c | < n x the dominant cell is tailed to the corner 
and compressing the buoyancy driven flow near the wall makes a tilt angle 9 ,, with respect to the vertical line E. 
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Increasing p c | , the LSC rotates around the vortex vector « max so that (),,, decreased and the buoyancy driven 

convection is expanded and the flow pattern between the wall and the edge, E is considerably simplified. 

Figure 4 shows that the intensity of the convection, \|/ max strongly depends on the rotation rate. Coupling to the 
radial velocity of the flow, the Coriolis force enhancement has a restrictive effect on the centrifugally forces, 
pronounced more at the lowest and the highest rotation rates (see figure 4), does not change the trend. 

Associated to the centrifugally driven flow, the vortex w max = (i/r){d V /Sr 2 + aV / 8z 2 -(Ur) dy/ /dr] 
magnitude is affected by the Coriolis and thermocapillary forces. 

The numerical results, with and without Marangoni flow, are plotted in figure 5. Increasing the rotation rate, 
co max is lowered to a minimum at around |n c |= 1.571 rad is (—15 rpm), and the raised sharply by further 

increase in F2 C . Again, thermocapillary forces modify the curve, particularly in the lower range of rotation, but 
the general feature of the curve remains the same. For the GaAs melt parameterized in the model by so high 

1/2 

Ma /PrGr as ~1 80, the result is unexpected. 

l/o 1/9 1/9 

According to the scaling analysis, in the interior of the melt y/ a ly/ p = 4 (Ma / Pr) /£ Gr , 
where y/ n and yr p are the magnitudes of Marangoni and buoyancy flow stream functions, respectively and c 
appears in at m = £ a 7' m . ix [29]. The temperature at the point Mj was computed for the applied range of rotation 
and found to be t m = 9.5 ± l .4 K leading to £ - 0.53 and yr a ly/ p ~ l/io . 

It was found that the tilt angle, 0 V is directly related to the Ekman shear angle in its model-adopted 
form X ( r ) = 4.746 x 10 _3 ^ 2 n c (r c 2 - r 2 ) , simply by tan 0 = 1 - tan X (r) at the point K (r = r c /2, z. — 8 E ) 

on the line E close to the crucible bottom. Herein, the angles are in degrees, n c in rpm, r c = 7.50 cm , and 
8 e = (v /|n t ,|) 1/2 is~1.0mmfor n c =-5.0 rpm as the Ekamn layer thickness. Increasing in n c , the number 

of spiral arms N E = | Re c I / 1 6 with separation &r = (i / 4N E )(r 2 / r) , increases. This can be readily verified 

that the line E separates apart the smaller spirals near the corner and the larger spirals near the axis of rotation, 
corresponding to the smaller and larger deformations of the flow volume element. Supposed that increasing the 
rotation rate intensifies the effect and referred to the simple relation between 0 lf/ and x (r) , development of the 

spiral shearing flow might be assumed to explain the LSC rotation around the vortex vector a max . 

3.3. Flow modes in the GaAs melt rotating with its container 

Referred to figure 3 (a-g), it appears that between the LSC and the rotation axis, there is an ovidal 
shape small roll with y/ /y/ max - 0.4 . The cell is located below the convex to melt part of the phase boundary 
and centered approximately on the line b crossing the inflection point I on the crystallization front. Emerged 
from the centrifugally driven cell, and circulating in the same direction as LSC, the intensity of this small 
convection found to be, in contrast to LSC, approximately independent of the rotation rate. Its shape, however, 
is modified by increasing £2 C . for lower rotation rates. In J< 7 .5 rpm and at the same shape, the thermal field 
exhibits an undulating structure in the melt leading to a radial non-uniformity of the thermal field close to the 
phase boundary. The effect is suppressed by further increase in the rotation rate (|n c |> 7.5 rpm ) due to the 

overstability known [20] as a distinctive feature of convection with rotation. This behavior implies that the 
Rayleigh-Benrad instability occurs at low rotation rates. In fact, when the buoyancy and Coriolis forces are 
comparable, (Gr/Ta ) - l , the Coriolis force will prevent vertical convection against the conservation of angular 
momentum and the oscillations decrease close to the onset boundaries. Herein, the ratio Gr/Ta equals to unity 
when |n c | = 8.78 rpm as confirmed by simulations of the fields in figure 3(a-f). 

For the imposed and calculated (at the point Mj) temperature differences AT r = a r m . a and 
0.8 < A T m j < 12 .5 K , respectively, the flow modes in the GaAs melt rotating at different rates are represented 
in figure 6 as the two-parameter regime diagram. Note that 8T Id z = 4.2 K / cm at the point M, found to be close 
to the data given by Kaoi et al [30], Therefore, the vertical temperature difference at is close to a T m i in the 
model. 
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The knee-shaped curve of transition obtained from the experimental results [15] for mercury ( 
Pr = 2.46 x 10 “ 2 ) shows that, at lower rotation rates the flow appears to be axisymmetric in the regime diagram. 
As mentioned before, the flow field is dominated by the Kiippers-Lortz instability for |n c |< 8.77 rpm . The 

enhancement of Coriolis forces complicates the flow pattern due to the baroclinic instability. Numerical results 
in figures 7(a) and 7(b) reveal that how the local temperature difference and velocity components at the mid- 
depth/mid -radius of the melt, namely at the point M, in the interior, vary with rotation. For more realistic case in 
the model, that is when Ma /0, both « m ; and T ln , decrease considerably with Q c in the lower range of 

rotation 5 < |q c | < 10 rpm . Remarkably, the flow velocity u m l is raised by a further increase in n c up to a 
maximum at n c = -12.5 rpm and then lowered more sharply in higher range of rotation, I2.5<|n c |<i5.0 rpm 
. Decreasing T m ; in the lower rotation range corresponds to higher value of the Rossby number Ro t (due to 

larger a T m t = T c - T m ; ) and lower values of the Taylor number ( Ta < 10 7 ) in the regime diagram. By further 
increase in n c , the onset of a baroclinic wave occurs at Ro T = 2.06 (at = 18 K and n c = - i .2 rpm ) 

corresponding to Ta = l .417 x 10 7 . For silicon melt, this critical Rossby number was 2.4 according to 
theoretical investigations [31]. The points which lie on the line II (see figure 6) stand for non-axisymmetric flow 
with regular wave [15], At rotation rates |n c |> 15 rpm and for the imposed at = A T max , the points located 

under the line II represent the stabilization of the flow at Ro T < 0.476 and Ta > 6.152 x 10 7 . The transition is 
qualitatively similar to experimental data [32] in which it occurs at Ro T number considerably lower than unity. 




Ta 



Figure 6: Flow modes at different rotation rates ( 2.5 <n c <20 rpm ) for the imposed horizontal and vertical 
temperature differences AT r = is k and 8.0 < at. < 12 .5 k , respectively. The transition curve obtained from 
Fein-Pfeffer’s experimental results [15]. The two lines stand for the upper and the lower boundaries for the 
occurrence of baroclinic instability. 
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Figure 7(a): Rotational effect on temperature variation in the interior of the melt. 




Figure 7(b): Rotational effect on the components of flow maximum velocity in the interior of the melt. 

IV. CONCLUSION 

The GaAs flow field structure and temperature distribution were described in the present model calculations. It 
was demonstrated that even perfectly axisymmetric geometries may exhibit flow instabilities. Simulations of the melt flow 
even at the lowest applied rotation rate fi t . = -0.524 rad Is revealed that buoyancy driven flow characterized by a 
unicellular meridional circulation is compressed by a strong centrifugally driven roll (LSC). Throughout the simulations, 
LSC appeared as the dominant flow in the multicell structure of the GaAs melt rotating with its container. In contrast to the 

1/2 

large magnitude of Ma /Pr Or , it was shown that the convective flow pattern does not depend strongly on the presence 
of thermocapillary forces. In fact, in the interior of the melt yr p ly/ a as the ratio between the buoyancy and 
thermocapillary forces largely exceeds the unity. 

Increasing the rotation rate 0.524 < |q c | < 2.094 rad Is , the dominant cell (LSC) found to rotate around the 
vortex vector ra max by an angle 9 ¥ directly related to the Ekman shear angle. Decreasing the intensity of the cell, the 

Coriolis force enhancement suppresses the undulating structure of thermal field. The so-called overstability of the melt was 
broken by further in the rotation rate. The appearance of the baroclinic wave region in the flow regime diagram corresponds 
to the transition from axisymmeteric to non-axisymmetric flow. The transition found to be sensitive to the applied 
temperature difference and as AT increases, a regular wave gives way to an irregular wave. The onset of baroclinic wave 
found to occur at Ro T = 2.07 when AT = 18 K . According to the experimental analysis [32] the transition is associated 
with significant increase in the amplitude of temperature fluctuation close to the phase boundary. As shown in the regime 
diagram, the baroclinic thermal wave disappears when the rotation-driven convection overcomes the buoyancy-driven 
convection. 

The present model analysis may ensure axisymmetric flow and thermal fields in the melt characterized 
by Ro r < 0.476 and Ta >6.152 x io 7 . The results correspond to the rotation rate |n c | > 1.571 rad is and 
// = |n c /n Jt |>3, close to the parameters used in the industrial scale growth ( |n c | = 1 .676 rad / s and ^ = 2.1 ) 
of the GaAs crystals. Notable is that the vortex of the dominant circulation <a max found to meet its minimum 
value at around |n c | = 1 .571 rad is. 
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